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Abstract 

We examine conditions for finite-time collapse of the solutions of the higher-order nonlinear Schrodinger 
(NLS) equation incorporating third-order dispersion, self-steepening, linear and nonlinear gain and loss, and 
Raman scattering; this is a system that appears in many physical contexts as a more realistic generalization 
of the integrable NLS. By using energy arguments, it is found that the collapse dynamics is chiefly controlled 
by the linear/nonlinear gain/loss strengths. We identify a critical value of the linear gain, separating the 
possible decay of solutions to the trivial zero-state, from collapse. The numerical simulations, performed for 
a wide class of initial data, are found to be in very good agreement with the analytical results, and reveal 
long-time stability properties of localized solutions. The role of the higher-order effects to the transient 
dynamics is also revealed in these simulations. 
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1. Introduction 


The present paper deals with the dynamical properties of a higher-order nonlinear Schrodinger (NLS) 
equation, expressed in dimensionless form: 

Ut + - i\ufu = 'yu + 5\ufu + Pu^xx + K\ufu)x + - iur) u{\u\‘^)x, (1.1) 

for all X G M, and t G [0,T] for some T > 0. In Eq. (O), u{xA) is the unknown complex field, subscripts 
denote partial derivatives, 7 , (5, cfr and /i, v are real constants, while 5 = 1 corresponds to normal group 
velocity dispersion. We supplement Eq. 0 with the periodic boundary conditions for u and its spatial 
derivatives up to the second order: 


u{x 2LA) = and 


d^u 

dxi 


{x + 2 L, t) 


d^u 

dxi 


j = l,2, 


Mx G K, Vi G [0,T], 


where L > 0 is given, and with the initial condition 


( 1 . 2 ) 


u{x^ 0 ) = U{){x)^ yX e . 


(1.3) 


also satisfying the periodicity conditions ( |l. 2 [ ). Note that the boundary conditions ( 1 . 2 ) will be incorporated 
in the Sobolev spaces of periodic functions, that will be used for the analysis of the problem 0 - 0 . 
Their definition will be recalled and discussed below. 

Equation 0 finds important applications in distinct physical and mathematical contexts, with the 
most prominent being that of nonlinear optics. Indeed, such models have been used to describe ultrashort 
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pulse propagation in nonlinear optical fibers in the framework of the slowly-varying envelope approximation 
niiais] (see also reviews [IJ [5] and references therein for recent work in few-cycle pulse propagation). In 
such a case, t and x in Eq. (0 play, respectively, the role of the propagation distance along the optical fiber 
and the retarded time (in a reference frame moving with the group velocity), while u{x,t) accounts for the 
(complex) electric field envelope. The propagation of optical pulses (of temporal width comparable to the 
wavelength) is accompanied by effects, such as higher-order (and, in particular, third-order) linear dispersion 
(characterized by the coefficient /3), as well as nonlinear dispersion (characterized by the coefficient /r). The 
latter is known as “self-steepening” term, due to its effect on the pulse envelope. These conservative higher- 
order terms (including the one of strength u) are naturally important not only in optics, but also in the 
contexts of nonlinear metamaterials liiTiii] and water waves laiioiiii]. All the above effects, which are 
not incorporated in the standard NLS equation {'j = S = P = jii = u = aR = 0), become important in the 
physically relevant case of ultrashort pulses. 

On the other hand, and in the same nonlinear optics context, for relatively long propagation distances, 
dissipative effects, such as linear loss (7 < 0) [or gain (7 > 0)], as well as stimulated Raman scattering (of 
strength > 0 ), are also included [Il[2j[3]. At the same time, it is also relevant to include the nonlinear 
gain {S > 0) [or loss {S < 0)] to counterbalance the effects from the linear loss/gain mechanisms [ 12 ]. For 
relevant physical settings, the role of such terms is to contribute to the possible stabilization of solitons. 
Notice that such solutions exist in the standard NLS equation, which possesses solitons of two different 
types: bright solitons (for s = —1) that decay to zero at infinity, and dark solitons (for s = + 1 ) that have 
the form of density dips on top of a continuous-wave (cw) background accompanied by a phase jump at 
the density minimum. Below we focus on the latter case of normal dispersion {s = +1), for which the 
linear/nonlinear dissipative effects usually lead to the decay of the cw background; this way, the composite 
object, composed by the background and the soliton, is destroyed. 

Dark solitons were first predicted to occur as early as in the early 70 ’s [ni[Ti], and since then have been 
studied extensively, both in theory and in experiments, mainly in nonlinear optics m and Bose-Einstein 
condensation o na n]. Concerning the higher-order NLS ( | 1 . 1 [ ), in Ref. [T9| it was shown, via a newly 
developed perturbation theory m, that dark and grey solitons can exist and feature a stable evolution 
-at least for certain parameter regimes and time intervals, within the validity of the perturbation theory. 
However, their long-time dynamics are under question. The results of Ref. m suggest that these dynamics, 
and different scenarios corresponding to various parameter regimes, can be very rich. 

The structure of the presentation and main findings of the present work are as follows. In Sec. we 
explore the finite-time collapse for the solutions of Eqs. (fnt-doi. Our analysis is motivated from the 
results and questions posed on the dynamics of the soliton background discussed in Ref. m- The discussion 
ther ein, shows that there can be a finite-time for which the background exhibits blow-up. Also in Section 
2 . 1 , based on energy arguments carried out on a spatially averaged L^-energy functional, we identify that 
the collapse behavior of the considered L^-functional is chiefly characterized by the gain/loss parameters 
7 and S. Regimes of collapse and decay of solutions are defined by 7, > 0 and 7, ^ < 0. However, for 

7 < 0, J > 0, we obtain the existence of a critical value on the linear gain. This critical value, say 7*, is 
separating finite-time collapse from the decay of solutions. An important conclusion regarding the regimes 
7, (5 > 0 and 7, (5 < 0 is related to the fate of localized solutions, with an initial profile resembling that of 
bright and dark solitons: referring to their long-time stability, such wave forms are unstable, and they may 
survive only for relevant short times. 

The above picture on the asymptotic behavior is verified completely by the results of the numerical 
simulations presented in Section We examine the evolution for three distinct types of initial conditions, 
including decaying and non-decaying data. In the case of cw initial conditions, we observe that the analytical 
upper bound is not only fulfilled, but also saturated by the analytical blow-up times. Furthermore, the 
numerical collapse time is independent of the strengths of the rest of the dissipative/conservative higher- 
order terms. The observed saturation is completely justified by the fact that the timescales for collapse of 
the background, and the analytical upper bound for the collapse of the functional, coincide in the case of 
the cw initial condition. Also, the evolution of the L^-functional is proved to be that of the cw-background. 
Concerning the critical value 7*, it is verified numerically that its analytical estimation is very accurate; for 
7 < 7* the solution decays and for 7 > 7* the solution blows-up in finite-time. In the case of the decaying 
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initial data (which are taken in the form of a sech function - as reminiscent of a bright soliton of the focusing 
NLS), we have observed very interesting transient dynamics prior to collapse when 7 > 7*, or decay when 
7 < 7*, characterized by a self-similar evolution of the initial pulse for small values of the parameters. As 
before, the analytical upper bound for the blow-up time is in agreement with the numerical one. However, 
it is not saturated since the evolution of the functional is dominated by the evolution of the amplitude 
of the localized solution. We have also tested numerically the behavior of the blow-up time with respect 
to large variations of the strengths of linear and nonlinear parameters P and cfr^ and found the following 
decreasing step function effect: the numerical blow-up time jumps to a smaller numerical value but remains 
constant within specific regimes of variation for the increasing parameters. The jump is associated with a 
transition to a different type of transient dynamics prior to collapse or decay. 

The numerical results were also found to be in excellent agreement with the analytical ones in the case 
of the non-decaying initial data -taken here in the form of a tanh^ function, which resembles the intensity 
profile of a dark soliton of the defocusing NLS, and is compliant with the periodic boundary conditions ( 1 . 2 ). 
The numerical blow-up time is in agreement with -and is also found to almost saturate to- its analytical 
upper bound, as in the case of the cw-initial condition, at least for small variations of the parameters, for 
reasons that will be discussed later in the text. In regimes of larger parameter values, we observe again 
the decreasing step function effect for the numerical blow-up time associated with a drastic change in the 
transient dynamics. In the case 7 < 0, (5 > 0, it appears that the condition j6 < 0 is not sufficient to define 
a stabilization regime with nontrivial dynamics or avoiding collapse. It will be shown in a subsequent work 
[ 21 ], that such a stabilization regime associated with emergent dynamics different from collapse or decay, 
may be defined only when 7 > 0 , < 0 . 


2. Conditions for collapse 


2.1. Collapse in finite-time for an L‘^-norm funetional 

The question of collapse concerns sufficiently smooth (weak) solutions of equation ([H}. The existence 
of such solutions, is guaranteed by the following local existence result associated with the initial-boundary 
value problem ([Tl]i-([Ot- The methods which are used in order to prove such local existence result in 
the Sobolev spaces of periodic functions Esum, are based on the lines of approach of [22]-[25]. The 
application of these methods to the problem (!.!)-(1.3), although involving lengthy computations, is now 
considered as standard. Thus, we refrain from showing the details here, and we just state it in: 


Theorem 2.1. Let uq G for any integer k>2, and /3, 7 , (5,/i, z/, (7^^ G M. Then there exists T > 0, 

sueh that the problem has a unique solution 

« G C([0, T], ^7(fi)) and G C([0, T], (^^)). 


Moreover, the solution u G depends eontinuously on the initial data uq G i.e., the solution 

operator 


S{t) : 

Uo 


^ te[0,T], 

1 -^ S{t)uo = u, 


( 2 . 1 ) 


is eontinuous. 

Here, denote the Sobolev spaces of 2L- periodic functions on the fundamental interval ft = 

[—L,L]. Let us recall for the shake of completeness their definition: 

1 j 

HperiS^) = C, wand — G j = l,2,...,k\ 

0^ u 

u{x), and -^—{x) for j = I, 2,..., k — 1, are 2I/-periodic}. (2.2) 
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Since our analytical energy arguments for examining collapse require sufficiently smoothness of local-in-time 
solutions, we shall implement Theorem | 2 . 1 | by assuming that /c = 3, at least. As it follows from the definition 
of the Sobolev spaces ( 2 . 2 ), this assumption means that the initial condition i 4 o(x), x G ^2, and its spatial 
(weak) derivatives, at least up to the 2nd-order, are 2I/-periodic. Then, it turns out from Theorem 


2.1 


that the unique, local-in-time solution u{x^t) of ( 1 . 1 ) satisfies the periodic boundary conditions ( 1 . 2 ) for 
t > 0, and is sufficiently (weakly) smooth. Such periodicity and smoothness properties of the local-in-time 
solutions are sufficient for our purposes (see Theorem 2.2 and Appendix A, below). 

Next, we adopt the method of deriving a differential inequality for the functional 


g- 27 t pL 

M{t) = J ^ \u{x,t)\^dx, 


(2.3) 


and then, showing that its solution diverges in finite-time under appropriate assumptions on its initial value 
at time t = 0 ; see [Ml EH EH EH EO] and references therein. Note that the choice of this functional is not 
arbitrary; in fact, it is a direct consequence of the conservation laws of the NLS model (see Appendix A). 
For a discrete counterpart of this argument applied in discrete Ginzburg-Landau-type equations, we refer to 
jST] . For applications of these types of arguments in the study of escape dynamics for Klein-Gordon chains, 
we refer to [32] . 

Theorem 2 . 2 . For uq G k > 3, let S{t)uo = u G C([0, T^ax), ^per(^)) local- in- time 


solution of the problem (1.1)- O, with [0, Tmax) be its maximal interval of existence. Assume that the 
parameter (5 > 0 and that the initial condition uq{x) is such that M(0) > 0. Then, Tmax 'Is finite, in the 
following cases: 


(0 


27 


1 + 


7 


for 7 7 ^ 0 , and 
(ii) Tmax ! 


6M (0) 
7 > —6M (0). 

25M ( 0 ) ’ ^ ^ 


0 . 


(2.4) 

(2.5) 

( 2 . 6 ) 


Proof: For any T < T^; 


since k > 3, due to the continuous embedding C{[0,T], 


C([ Q, T], T^(fl)) [26] . the solution S{t)uo = u G C([0, T], T^(fl)). Furthermore, it follows from Theorem 


2.1 


that Ut G C([0, T], T^(0)). Then, by differentiating M(t) with respect to time, we find that 


dM{t) 

dt 




= -7- 


/: 


\u\‘^dx ■ 


o- 27 t 


-Re 


utudx. 


(2.7) 


In the second term on the right-hand side of (2.7), we substitute Ut by the right-hand side of Eq. 
Then, after some computations (see details in Appendix A), Eq. (2.7) results in: 


dM(t) 


= 5 - 




dt L 

Next, by the Gauchy-Schwarz inequality, we have 


/: 


\u\‘^dx. 


( 2 . 8 ) 


/: 


^dx < 


if- 


wy^dx 


(2.9) 


Therefore, for the functional M(t) defined in (2.3), we get the inequality 

.,-274 f rL \ 3 


M{t) < 




lul'^dx 


( 2 . 10 ) 
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which in turns, implies the estimate 


M{tf < 



for all t € [0,Tmax)- On the other hand, from (2.8) we have that 





L dM{t) 
S dt 


( 2 . 11 ) 


and hence, we may rewrite ( 2 . 11 ) as 




2 dM (t) 


26 dt 


or 


dM{t) 

dt 


2 — 


> 2Se^'>K 


Since I 


M(t) 


dM(t) 


, we get from (| 2 . 12 |) the differential inequality 


1 


d 


dt 


< -2Se^^\ 


Integration of (2.13) with respect to time, implies that 

1 1 

TTTT < T:rwT - 2<5 / 
M{t)- M (0) Jo 

and since M (t) > 0 , we see that M( 0 ) >0 satisfies the inequality 

26 f e^'^^ds < ^ 

Jo 


M(0)’ 


From (2.14), we shall distinguish between the following cases for the damping parameter 7 : 


( 2 . 12 ) 


(2.13) 


(2.14) 


• We assume that the damping parameter 7 7 ^ 0. In this case, (2.14) implies that 

^ (e^'yt - 1 ) < ^ or < 1 + _ 

27 ^ - M ( 0 ) ’ - (5M ( 0 ) ■ 

Thus, assuming that > —1, we derive that the maximal time of existence is finite, since 


t < 7 - log 
27 


1 


6M ( 0 ) 


The inequality above, proves the estimate of the collapse time (2.4) under assumption (2.5), that is, 
case (i) of the Theorem. 


Assume that the damping parameter is 7 = 0. Then, Eq. (2.14) implies that 26t < 


t < 


1 


26M (0)' 


This inequality proves the estimate of the collapse time (2.6) in the absence of damping, that is, case 
(ii) of the Theorem. □ 
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From condition (2.5) on the definition of the analytical upper bound of the blow-up time 


J7,(5,M(0)] = ;^log 
27 


1 


7 


6M (0) 


given in (2.4), we define a critical value of the linear gain/loss parameter as 

7* = -6M{0). 

We observe that 


lim = +00, 

while Tniax[7, ^(0)] is finite if 


(2.15) 


(2.16) 

(2.17) 


. 5 }< 

7 > 7 


(2.18) 


according to the condition (2.5). Then, (2.17) suggests that when > 0, the critical value 7* may act as a 
critical point separating the two dynamical behaviors: blow-up in finite-time for 7 > 7* and global existence 
for 7 < 7*. We shall investigate this scenario numerically. 


We also remark that the analytical upper bound for the blow-up time (2.6) in the case 7 = 0, 

T^e^[5,M{0)] ^ 


26M (0) ’ 


(2.19) 


is actually the limit of the analytical upper bound (2.15) for 7 > 0 as 7 ^ 0, e.g.. 


lim Tmax[75 ^5 ^(0)] — Tjnax [<5,M(0)]. 

7^0 


( 2 . 20 ) 


3. Conditions for collapse in finite-time: numerical results 

In this section, we verify numerically the collapse of different types of initial data. We shall study three 
particular cases for the initial data: (a) cw initial conditions, (b) sech-profiled initial conditions, resembling 
a “bright soliton”, and (c) tanh^-profiled initial conditions, resembling a “dark soliton”. 

All numerical experiments have been performed using a high accuracy pseudo-spectral method, incorpo¬ 
rating an embedded error estimator that makes it possible to efficiently determine appropriate step sizes. 
We denote as time of collapse the time where the numerical scheme detects a singularity in the dynamics. 
The numerical singularity is detected when the time step is becoming of order < 10“^^. Concerning the 
implementation of the periodic boundary conditions for cases (b) and (c) of the initial data, let us also 
recall that in both cases, these conditions are strictly satisfied only asymptotically, as L ^ 00; for a finite 
length L, these initial profiles, as well as their spatial derivatives, have jumps across the end points of the 
interval [—L, L]. However, these jumps have negligible effects in the observed dynamics, as they are of order 
exp(—L). More precisely, since the smallest value for L used herein is L = 20, the derivatives are 0(10“^), 
or less. 


3.1. Continuous wave initial conditions. 

The first example for our numerical study, concerns cw initial data, of the form: 

• Kttx 

uo{x) = ^ , K > 0, 


(3.1) 


of amplitude e > 0 and wa ve-nu mber K > 0. This type of plane wave dynamics can be understood by means 
of the substitution in Eq. (1.1) for ansatz of the form u{x^t) = W{t)e^^^^ in line with the above presented 
numerical simulations. Then the resulting dynamics is: 


Wt 


.sK‘^ 


W - i|iy|2iy = -fW + S\WfW - ipK^W + ifiKW 


(3.2) 
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Subsequent introduction of a phase factor in the form of W{t) = with uj = sK‘^12 — + jiK 

absorbs through this gauge transformation a number of the terms in Eq. (|3.2|) finally leading to: 


iw — \w\‘^w = i^w + iS\w\‘^w^ 

where w = dw/dt. Then, employing the polar decomposition w = we get 

h = yh + 5h^. 


(3.3) 


Assuming that the initial height of the wave background is h(0) = ho, the ODE (3.3) has the general solution 


h\t) = 


7 + hhg - • 


The solution of (3.4) blows-up in the finite-time 


= 7^ log 

27 


5hl, 


(3.4) 


(3.5) 


exactly under the same assumptions on the strengths 7 and 5 stated in Theorem 2.2 


Recalling the above analysis, for the initial data (3.1) we have 

= ^ r e^\e-^^\^dx = e\ 




Theorem 2.2 and (2.15) assert that the analytical upper bound of the collapse time 6, M{0)] is 


rmax[7><^>e] = Tlog 
27 




According to (2.16), the critical value 7* in this case [see also (3.6)] is: 

7* = 


Also, (2.19) asserts that the analytical upper bound for the collapse time in the case 7 = 0 is 




1 


2he2' 


(3.6) 


(3.7) 


(3.8) 


The left panel of Eig. illustrates the comparison of the analytical upper bound on the collapse time 


(3.6), against its numerical value, as a function of the parameter 7, fixing the rest of parameters at specific 


examples of values. We observe that the analytical upper bound seems to be a very accurate estimate of the 
numerical collapse time, and is almost saturated. This, one would argue, is not surprising provided that the 
solution preserves its plane wave structure, since the latter indeed saturates the Cauchy-Schwarz inequality 


(2.9) used to derive our differential inequality from which the collapse time was obtained. Indicative of 


the accuracy of the analytical estimate is the comparison between the analytical upper bound for the case 


7 = 0, and the numerical blow-up time for the considered set of parameters: the analytical value (3.8) 
is Tmax[0.1,1]=5 (marked by the blue square in the figure), and the numerical blow-up time (the red dot 


inside the square) is Tnum =4.99. The analytical upper bound ( |3.6[ ) does not depend on the wavenumber K 
and the half-length L of the symmetric spatial interval [—L,!/]. This independence is also preserved by the 
numerical blow-up time as illustrated in the right panel of Eig. The observed independence is justified by 
the same reasons as above (i.e., that the inequality of ( |2.9[ ) is indeed an equality for this type of solutions, 
provided they maintain their spatial form). 

Eigurej^ illustrates the comparison of the analytical upper bound on the collapse time (3.6) against the 
numerical collapse time, as a function of the parameter > 0, the amplitude of the plane wave e, and the 
parameters <7^^, /d, /i, and u. In the upper panel, we observe again the high accuracy of the analytical upper 
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Figure 1: (C olor Online) Left panel: Numerical blow-up time for the plane wave initial data {33 against its analytical upper 
bound ( |3.6| >, as a function of 7 G [—0.09,0.1], for 6 = 0.1, [i = v = 0.01, (Jr = 0.02, ^ = 1. Continuous (blue) curve depicts 
the analytical upper bound and (red) dots the numerical blow-up time. Here L = 100, while iF = 10, e = 1 for the cw. The 
analytical upper bound for 7 = 0 is marked by the (blue) square. Right panel: Independence of the blow-up time on the wave 
number K for cw initial data and on the length L of the spatial interval (inset). Here, 7 = —0.08, the rest of parameters are 
as above, and for fixed L = 100, K G [10,100] (main panel), while for fixed K = 10, L G [100,1000] (inset). For this set of 
parameters, the analytical upper bound for the blow-up time is Tmax = 10.05 (lines-marked by the black square in the upper 
figure), while the numerical blow-up time is Tnum ~ 10.01 (dots). 


bound for the collapse time ( 3 . 6 ), as a function of the nonlinear gain/loss parameter and as a function of 
the amplitude e. In the bottom panel, we illustrate the independence of the numerical blow-up time from 
the parameters / 3 , /i and z^, as predicted by the analytical estimate. We remark that the last two terms 
of the nonlinearity of <0) can been written as 


F[u\ = [{jj, + 7 - + [il^ + v) - iaR]u^Ux 

Thus, we have integrated the system varying the sum jj. -\- u. 

In Fig. 1 ^ we present the results of a numerical study testing the behavior of the critical value ( 3 . 7 ) 
for the linear gain/loss parameter 7, as a separatrix between global existence and collapse. The left panel 
shows the evolution of the density for plane wave initial data ( 3 . 1 ), for the parameters example of the upper 


panel of Fig. and the choice of the parameter 7 exactly at the critical value, 7 = 7* = — 0 . 1 . We observe 
the blow-up of the density, found numerically at time Tnum = 27 . 56 . The collapse is similar to that of an 
unstable cubic ordinary differential equation (ODE), with a rapid, spatially uniform, monotonic increase 
of the density as the evolution approaches the collapse time. Notice the preservation (both in this case 
and in the case discussed below) of the plane wave character of the solution. On the other hand, slightly 
perturbing the parameter 7 < 7*, we observe the complete change in the asymptotic behavior, exhibiting 
the decay of solutions. This is shown in the right panel, illustrating the evolution of the density for plane 
wave initial data for the parameters example of the upper panel, but for 7 — 0.101 < 7* = —0.1, slightly 
below from the critical value 7*. The decay is again spatially uniform, resembling this time, that of a stable 


cubic ODE. We observe that the time-scales for blow-up ( 2 . 15 ) —and ( 3 . 6 ) i n th e case of the plane wave 
initial data— coin cide with the time-scales for blow-up of the background h‘^ ( 3 . 5 ). Eurthermore, it follows 
from the solution ( 3 . 4 ) that lim^^oo ^(^)^ = 0 , if 7 < for J > 0 , i.e., the condition for decay of the 


background h{t) coincides with the condition 7 < 7* for the decay of the solution in the regime 7 < 0, > 0. 

This coincidence explains the saturation of the analytical upper bound ( 3 . 6 ) by the numerical blow-up times 


from a complementary perspective (to that of the saturation of the inequality of Eq. ( 2 . 9 )) and the overall 


agreement of the parametric dependencies between the numerics and the theoretical analysis above. 

3 . 2 . Decaying (sech-profiled) initial conditions. 

The second example for our numerical study concerns sech-profiled initial data 


uo{x) = esechx. 


( 3 . 9 ) 
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Figure 2^ _ (Color Online) Top left panel: Numerical blow-up time for the cw initial data against its analytical upper 

bound ( |3.6| , as a function of (5 G [0.001, lO j, fo r 7 = 0.04, /j, = 1 / = 0.01, aji = 0.02, (3 = 1. The continuous (purple) curve 
corresponds to the analytical upper bound j3.6| > and (red) dots to the numerical blow-up time. Here, L = 100, while K = 10, 
e = 0.1 for the cw. Top right panel: Numerical blow-up time (dots) against the analytical upper bound (continuous (black) 
curve) as a function of e G [0.001,1], for 7 = 0.04, 6 = [j, = u = 0.01, aji = 0.02, ^ = 1, L = 100, K = 10. Bottom left panel: 

Independence of the blow-up time on the parameters aji and /3 (inset). Here, 7 = 0.04, S = 0.3, e = 0.1, K = 10, L = 100; the 
analytical upper bound for the blow-up time (black circles) is Tmax = 33.28. (lines), is marked by the black circle in the top 
panels. In the main panel, (3 = 1, and an G [0.01,0.5], while aR = 0.3, /3 G [0.1,1.3] in the inset. The numerical blow-up time 
is Tnum ~ 33.21 (dots). Bottom right panel: Independence of the blow-up time on parameters /r + z/ G [0.01,0.5], for (3 = 1, 
(Jr = 0.3; the rest of parameters are as in the upper right panel. 



Figure 3: (Color Online) Left Panel: Collapse of the density |np for plane-wave initial data ( |3.l| , as found for d = 0.1, e = 1, 
and the critical value for 7 = 7*— 0.1. The rest of parameters are K = 10, L = 200, ^ = 0.01, u = 0.01, a = 0.02, ^ = 1. 
Right Panel: Decay of the density |np for plane-wave initial data, as found for 5 = 0.1, e = 1, and 7 = —0.101 < 7 * — 0.1. The 
rest of parameters are the same as in the upper panel. 
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Figure 4: (Co lor Online) Left panel: Numerical blow-up time for the sech-profiled initial data ( |3.9| > against its analytical upper 
bound ( |3.10| >, as a function of 7 G [0.1, 0.5]. Parameters S = 0.01, fi = u = 0.01, cfr = 0 . 02 , (3 = 0 . 02 . The dashed (blue) 
curve is traced for the analytical upper bound ( |3.10| > and (black) dots for the numerical blow-up time; here L = 50. The 
inset illustrates the numerical blow-up time (dots) against the analytical upper bound (dashed (blue) curve) as a function of 
(Jr G [0.4,1]. Parameters 7 = 0.1, 5 = 0.01, = u = 0.01, (Tr = 0.02, j3 = 0.02, L = 50. Right panel: Numerical blow-up 

time for the sech-profiled initial data ( |3.9| > against its analytical upper bound ( |3.10| >, as a function of L G [20,140]. Parameters 
7 = 0.01, 5 = 0.01, = u = 0.01, Jr = 0.02, [3 = 0.02. The dashed (red) curve is traced for the analytical upper bound ( |3.10| > 

and (black) dots for the numerical blow-up time. The insets illustrate the numerical blow-up time (dots) against the analytical 
upper bound (dashed (blue) curve) as a function of /3 G [0.01, 0.2] (left inset) and /3 G [0.4,1] (right inset) . Parameters 7 = 0.1, 
5 = 0.01, fi = u = 0.01, (JR = 0.02, (3 = 0.02, L = 50. 


of amplitude e > 0. Such iuitial data correspoud to the profile of a “bright solitou” as au iuitial state. lu 
this case, 


M(0) = i/ 


^2 1 2 7 e^tauhL 

e sech xdx = ---, 

L ^ 


aud Theorem 2.2 aud (2.15) assert that the aualytical upper bouud of the collapse time M{0)] is 


^max['T: ^5 ^5-^] — q 
2-f 


1 + 


7 T 


(^e^tauhl/ 


Also, (2.19) asserts that the aualytical upper bouud for the collapse time iu the case 7 = 0 is 


L 


T^^[S,e,L] - 2j^2tanhi:' 


Accordiug to (2.16), the critical value 7 * iu this case [cf. (3.10)] is: 

^ Se^tauhL 


L 


Moreover, it follows from (3.12), that for a sufhcieutly large half-leugth parameter L, 


7 


0 . 


(3.10) 


(3.11) 


(3.12) 


(3.13) 


Iu the left pauel of Fig. we preseut the uumerical study ou the comparisou of the uumerical blow-up 
time for the sech-profiled iuitial data (3.9), agaiust its aualytical upper bouud (3.10), as a fuuctiou of 
7 G [0.1, 0.5]. The rest of parameters are 6 = 0.01, fi = u = 0.01, cfr = 0.02, jS = 0.02 aud L = 50. We 
observe that the aualytical upper bouud is agaiu satisfied. However, it is uot saturated as iu the case of 
the cw iuitial data (as discussed above). Here, the time of collapse differs systematically by about 15-30%, 
suggestiug the poteutial for more iuterestiug dyuamical evolutious that we uow explore. Iu fact, the case of 
the sech-profiled iuitial data, is much differeut siuce the iuitial data are spatially localized. By coutiuuity, 
the solutiou will be also spatially localized for a fiuite iuterval [0, r] prior to collapse. Iu this case, for the 
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Figure 5: (Color Online) Dynamics for sech-type initial data, for /3 = 0 . 02 , S = 0.01, e = 1, jj, = 0.01, ly = 0.01, a = 0 . 02 , 
and L = 50. Upper left panel: Evolution of the density |np towards collapse, when 7 = 0.1 > 7 * ~ 0. Upper right panel: 
Decay of the density |np when 7 = —0.1 < 7 * 0. Bottom panel: Evolution of |n(0,t)| for the center of the pulse towards 

collapse when 7 = 0.1 > 7 * 0, corresponding to the density evolution shown in the upper left panel. The insets demonstrate 

snapshots of the density |np for specific values of |n( 0 , t)|-marked by lines connecting the density profiles. 


height of the wave-background we have h(t)^ = 0 = \u{x,t)\‘^ —the minimum of the density of 

the solution. Also, in this case we have the inequality 







7 

< 

M(t) 

^ j 



g-27t 



< 

2L 

max \u 
xe[-L,L] 

{x, f)|^ 21/ 


= e max 

\u{x,t)\‘^, 



xe[-L,L] 


(3.14) 


satisfied by the amplitude ^ \u{x,t)\‘^ of the localized solution, for all t G [0,t]. Thus, in the case 

of spatially localized data, we expect that the collapse will be predominately manifested as the increase of 
the amplitude, which grows faster than the averaged norm M(t), due to (3.14). Consequently, the amplitude 
may be expected to blow-up at an earlier time than M{t). 

It is also interesting to remark the slow logarithmic increase of the numerical blow-up time with respect 
to the parameter L G [20,140], as predicted by the analytical upper bound (3.10), and shown in the right 
panel of Fig. |^i.e., although the quantitative collapse time may be different, qualitatively our expression 
very accurately captures the relevant functional dependences. The parameter values for this example are 
7 = 0.01, 5 = 0.01, ^ = u = 0.01, gr = 0.02, f3 = 0.02. 

The insets in both panels of Fig. [^examine numerically the independence of the collapse time on others 
among the rest of the nonlinear parameters. In both examples we have considered the analytical upper 
bound Tmax = 31, for 7 = 0.1, 6 = 0.01, fi = u = 0.01, gr = 0.02, P = 0.02, and L = 50. The inset in the 
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Figure 6 : (Color Online) Dynamics for sech-type initial data, when /3 = 0.2. The rest of parameters are fixed as in Fig. 
Upper left panel: Evolution of the density |np towards collapse, when 7 = 0.1 > 7 * ~ 0. Upper right panel: Decay of the 
density |np when 7 = —0.1 < 7 * 0. Bottom panel: Snapshots of the density |np towards collapse when 7 = 0.1 > 7 * 0, 

corresponding to the density evolution shown in the upper left panel. 


left panel shows the small variations of the numerical blow-up time with respect to cfr G [0.4,1], around 
a mean constant value Tnum,crR ~ 24. The left inset in the right panel shows a ’’decreasing step function” 
effect for the numerical blow-up time with respect to P G [0.01,0.2], while the right inset shows that the 
numerical blow-up time has reached a mean constant value Tnum,/3 ~ 17.8 when (3 G [0.4,1]. The difference 
between Tnum,aR and Tnum,/3 suggests that the real blow-up time may jump to different numerical values as 
the nonlinear parameters vary (the behavior for /i, is similar), although remaining constant within 

specific intervals of variation. This effect may be related to transition between different dynamics of the 
initial pulse prior to collapse, as will be shown below. 

In Fig.|^ we present the results of a numerical study verifying that (3.13) is an approximate critical point 
separating global existence from collapse. Parameter values used for this study are P = 0.02, S = 0.01, e = 1, 
ju = R = 0.01, cfr = 0.2 and L = 50. The left upper panel illustrates the evolution of the initial density ( |3.9[ ) 
for t G [0,15], with 7 = 0.1 > 7* « 0, towards collapse at time Tnum = 27.85 (< Tmax = 31 —the analytical 
upper bound). We observe a nontrivial collapse mechanism, characterized (after a transient stage) by a 
self-similar evolution. Initially, the pulse amplitude is decreased but its width is increased (i.e., the pulse 
tends to disperse); nevertheless, as the gain mechanisms start playing a crucial role in the dynamics, the 
effect on the amplitude is reversed, and both width and amplitude grow self-similarly towards a collapsing 
state. To highlight this behavior, we show in the bottom panel the evolution of |iz(0,t)| of the center of 
pulse for t G [0,20]. The amplitude decays and after reaching a minimum, it starts growing. The insets 
are showing specific snapshots of the density \up corresponding to specific times and values of norm. The 
first one, shows the initial condition at t = 0, while the second shows the density profile at t = 4, where the 
amplitude has decreased and the width of the pulse has increased. The last one shows the density profile 
at t = 13, where both amplitude and width have increased, and the pulse evolves towards collapse. On the 
other hand, as shown in the upper right panel of Fig.|^ when 7 = —0.1 < 7* the initial density decays at an 
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exponential rate. This decay still occurs in a self-similar manner, where the amplitude of the pulse decays 
faster than the increase of its width. Similar dynamics towards collapse (or decay) have been observed for 
small variations of /3, as well as of the rest of the nonlinear parameters. 

Larger variations of these parameters may lead to a different type of transient dynamics towards collapse 
or decay. In Fig. we present the study for (3 = 0.2 and the rest of parameters fixed as in Fig. The 
top left panel shows the evolution of the initial density (3.9) for t G [0,12] towards collapse. The numerical 
blow-up time is decreased further to Tnum = 20.89, and the dynamics is very different from the self-similar 
evolution of the previous case, due to the considerable dispersion effect. This difference is highlighted in the 
bottom panel, showing snapshots of the density profiles at specific times, prior to collapse. The initial pulse 
is decomposed to a “core pulse” of decreased amplitude and a nonlinear wave train, both traveling to the 
right at slow speed (upper row). Progressively (bottom rows), and as the gain becomes more prominent, 
the amplitudes of the core and of the decomposed waves start increasing (note the change of scaling in the 
amplitude axis). While at a particular time instant, the amplitude of the wave-train is almost reaching that 
of the core, as shown in the snapshot for t = 10.5, eventually the core gains amplitude in comparison with 
the rest of the waves, as shown in the snapshot for t = 11. A sudden increase of the amplitude of the core 
is observed at t = 12, leading to collapse. For further variations of the parameter /3, the numerical blow-up 
time is stabilized —an effect illustrated in the inset of the right panel of Fig. EF and the collapse dynamics 
were found to be similar. The top right panel present the results of the numerical study for 7 = —0.1 < 7*. 
The initial pulse is again decomposed to a wave train, but now the dynamics are reversed to decay instead 
of collapse. 


3.3. Non decaying (tanh? -profiled) initial conditions. 

The third example of initial condition which we shall examine is the case of tanh^-profiled initial data. 


uo{x) = etanh^x. 


(3.15) 


of amplitude e > 0 (here, the discussion on the values of the derivatives on the boundaries should be 
recalled). Such initial data resemble a density (|r^P) profile reminiscent of the intensity dip of a dark soliton 
of the defocusing NLS, on the top of the cw background, which is compliant with our periodic boundary 
conditions. Note, that a simpler tanh-type initial condition (real or complex to model the black or gray 
solitons of the NLS) would violate our periodic boundary conditions (and would require a different type of 
boundary conditions, such as Neumann). While the observations here are still valid for the latter case, we 
will report corresponding results in a future communication. For the initial data ( |3.15 ), 


M(0) = i / 


L 

/ e^tanh^xdx = 


31/ — (3 + tanh^l/)tanhl/ 
3L 


Thus, in the present case. Theorem 2.2 and (2.15) assert that the analytical upper bound of the collapse 
time Tmax[7:^,^(0)] is 


= 7^ log 
27 


1 


37I/ 


6e‘^[{3L — (3 + tanml/)tanhl/] 


Then, (2.19) asserts that the analytical upper bound for the collapse time in the case 7 = 0 is 

T^^[5,e,L]= 


25e^[{2)L — (3 + tanh^L)tanhL] 


Furthermore, Eq. (2.16), implies that the critical value 7* in the present case [see (3.16)] is: 

3i — (3 + tanh^l/)tanhl/ 


7* = 


3L 


(3.16) 


(3.17) 


(3.18) 
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It is interesting to recover that for a sufficiently large L, the initial value of the functional M(t) at t = 0 


(3.15), the analytical upper bounds for the collapse time (3.16) and (3.17), and the critical value 7 * given 
in (3.18), approximate the relevant values for the example of the plane-wave initial conditions discussed in 


Sec. 3.1. Namely, for large L, 


M(0) 


j7,(5,e,L] « ;^log 
27 


1 + 


2 _ 


for 7 7 ^ 0 , 




25e2’ 


for 7 7 ^ 0 , and 7 * ^ —8e^ 


(3.19) 

(3.20) 


This result is of course not surprising given that the fact that as L increases the weight nature of the average 
in the definition of M( 0 ) yields a progressively vanishing contribution from the localized part of the relevant 
initial condition profile. The approximations (3.19) and (3.20) may indicate that the numerical collapse 


times may not only respect, but also be much closer to their analytical upper bounds as in the case of the 
plane-wave initial data. Another feature, interesting by itself, which may support the above argument, is 
that the collapse may be manifested by the collapse of the wave background for small values of the dispersion 
parameter P and of the nonlinear parameters, as discussed in Sec. 3.1. This fact should be also supported by 
the expectation that the dip of the density profile should have a little effect to the averaged norm M(t), and 
the main contribution should come from the wave background. More precisely, the approximate equation 
on the evolution of M{t) should be also valid in the case of the initial data (3.15). 
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Figure 7: (Color Online) Upper panel: Snapshots of the evolution of the density |np towards collapse for tanh^-type initial 
data, when 5 = 0.02, e = 1 and 7 = 0.025 > 7 * = —0.02. The rest of parameters are L = 50, /r = 0.01, u = 0.01, a = 0.02, 
/3 = 0.01. Bottom left panel: An upper view for the collapse of the density |np for tanh^-typ e initi al data for the set of 
parameters of the upper panel. Bottom right panel: Comparison of the analytical upper bound ( |3.16| > (dashed (red) curve) 
against the numerical blow-up time (dots) as a function of 7 . Rest of parameters are fixed as in the study of the upper panel. 


In the top and bottom left panel of Fig. we present the numerical results on the evolution of the 
density of initial data (3.15), for e = 1, 7 = 0.025 and 5 = 0.02. The rest of parameters are /i = 0.01, 
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V = 0.01, G = 0.02, ^ = 0.01, and the system is integrated for L = 50. For the given values of parameters, 
the critical value 7 ^ —0.02 < 0, and since 7 > 7 *, we are in the regime of finite-time collapse. The top 
panel shows snapshots of the density profiles at specific times, and the bottom left panel shows a (bottom) 
view of the density evolution, prior to collapse. The increase of the density is manifested by the increase of 
the height of the wave-background and accompanied by the nontrivial dynamics of the initial dip, prior to 
collapse at numerical time Tnum = 15.9, while the analytical upper bound is T^ax = 16.5. As highlighted in 
the upper panel, the initial density dip appears to separate into two traveling (at low speed) dark pulses of 
progressively increasing depths, emitting radiation in the form of linear waves of increasing amplitude. In 
fact, the emerging patterns bearing an oscillatory wavefront, together with a solitary wave are reminiscent 
of dispersive shock waves (DSWs) [33]. The generation of two dark soliton pulses can be understood by the 
fact that the initial condition has a trivial phase profile (i.e., does not have the phase jump needed for a 
dark soliton); thus, the initial dip separates into two dark solitons moving to opposite directions so that the 
phase of one soliton is compensated by the other soliton, with the total phase retaining its initial asymptotic 
structure - see, e.g.. Ref. |34|. On the other hand, the clearly observed asymmetry in the evolution profile 
is apparently produced by parity breaking terms, such as the first and third spatial derivatives in Eq. 

The bottom right panel of Fig. shows the comparison of the analytical upper bound (3.16) (dashed (red) 
curve) against the numerical blow-up time (dots) as a function of 7 , when the rest of parameters are fixed 
as in the study of the upper panel. We observe that the numerical blow-up are quite close to their analytical 
upper bounds, similarly (for the reasons indicated above) to the plane wave case. 



Figure 8 : (Color Online) Left panel: An upper view for the decay of the density |np for tanh^-type initial data, as found for 
S = 0 . 02 , e = 1, and 7 = —0.025 < 7 * = — 0 . 02 . The rest of parameters are L = 200, /r = 0.01, u = 0.01, a = 0 . 02 , /5 = 1. 
Right panel: A bottom view for the decay of the density |np for tanh^-type initial data. Parameters as in the upper panel. 


In Fig. keeping the rest of parameters as above and integrating the system for L = 200, we present 
the results of the numerical study for slightly perturbing the parameter 7 from its critical value 7 * = — 0 . 02 , 
falling within the decay regime, e.g, by considering 7 = —0.025. The left (right) panel shows a top (bottom) 
view of the evolution of the density manifested by the decay of the height of the wave-background again 
accompanied by nontrivial dynamics of the initial dip, prior to decay: it separates again into two traveling 
(at low speed) dark pulses, this time of slowly decreasing depths. The two decaying pulses are emitting 
radiation again, which also decays, following the decay of the wave-background. Once again, on both sides 
the oscillatory waveforms coupled to the dark solitary waves are somewhat reminiscent of DSW patterns. 

For larger valued-parametric regimes of the higher-order effects, the transient dynamics prior to collapse 
or decay were found to be potentially quite different. Concerning collapse, we have recovered the same 
effect for the variation of the numerical blow-up time associated with the transition to a different type of 
collapsing dynamics. In Fig.|^ we present snapshots of the density evolution for gr = 0.2 and the rest of the 
parameters fixed as in the case of Fig. —the case of collapse. While the initial dark pulse is still separated 
into two dips as explained above, we also observe stronger “one-sided” radiation to the left, due to the 
third-order dispersion, which is increased due to the (relatively strong) Raman effect. The two dips and the 
right-part of the background seem to have a similar dynamics as in the case of Fig.[^ however the collapse is 
manifested by the sudden increase of the amplitudes of the oscillatory waveforms in the background formed 
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to the left of the dips. 
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Figure 9: (Color Online) Snapshots of the evolution of the density |np towards collapse when 7 = 0.025 > 7 * = — 0 . 02 , all the 
parameters fixed as in the study of Fig. but for the increased value of aji = 0.2. 


The transient dynamics of our system can be very rich. As an example, in Fig. we show a contour 
plot of the space-time evolution of the density \u\‘^ in the decay regime 7 = —0.025 < 7* = —0.02 (the 
rest of parameters fixed as in the study of Fig. [^, but for increased dispersion effect, namely ^ = 0.2. 
The initial dip is decomposed to an almost spatially periodic pattern, traveling to the left. These periodic 
structures, in turn, have different survival times prior to the their final decay. A more detailed study of such 
phenomenology is warranted although it is outside the scope of the present investigation. 


4. Conclusions 


In this work, analytical studies corroborated by numerical simulations, considered the collapse dynamics 
of the higher-order NFS equation (1.1) supplemented by periodic boundary conditions. The model has 


important applications in describing the evolution of ultrashort pulses in nonlinear media (such as optical 
fibers, metamaterials, water waves), where linear and nonlinear gain/losses and higher order conservative 
and dissipative effects are of importance. We have shown by means of simple analytical energy arguments 
and verified numerically, that the collapse/decay dynamics are chiefly dominated by the linear and nonlinear 
gain/loss of strength 7 and S respectively. We identified a critical value 7* in the regime 7 < 0, (5 > 0, 
which is separating finite-time collapse from the decay of solutions. An important conclusion regarding this 
regime, is related to the fate of localized solutions such as pulse and dip-type waveforms [19]: referring to 
their true long-time dynamics, such states only survive for relatively short times, as they are subject to 
either catastrophic collapse or to eventual decay for generic parameter values. 

The results of the numerical simulations for a variety of initial conditions suggest a very good agreement 
between the analytical estimates for the collapse times and the numerical ones, as well on the accuracy of 
the critical value 7*. This agreement was excellent for plane wave initial data, and very good for waveforms 
bearing a density dip, while it was less good in the case of initial data resembling bright solitary wave 
structures, for reasons explained in the text. Furthermore, the numerical simulations gave insight on the 
role of the strengths of the higher order effects on the transient dynamics prior to collapse or decay, which 
were found to be essentially different. They have also indicated that these strengths may seriously affect 
the survival times of the solutions. 

Our results pave the way for a better understanding of the dynamics of the considered higher-order NFS 
equation. An interesting direction is to discuss the collapse dynamics for the case of s = — 1. It is also 
certainly relevant, to explore the associated dynamics in the stabilization regime 7 > 0, (5 < 0 respectively 
(see [21] for some recent results), depending on the higher-order effects and the type of initial data. Future 
investigations should also consider the possibility of collapse of higher-order norms, and the effect on this 
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Figure 10: (Color Online) Contour plot of the evolution of the density |np for 7 = —0.025 < 7 * and /3 = 0 . 2 . The rest of the 
parameters is fixed as in the the example of Fig. 


type of collapse of the higher-order nonlinear terms, as well as of higher-dimensional settings, which may 
further modify the dynamical regimes analyzed herein. 
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Appendix A. Derivation of Eq. (2.8) 


Here, for the sake of completeness, we provide details for the derivation of Eq. (2.8) and also explain the 
choice of the functional M{t) [cf. Eq. (2.3)]. 


First we note that the second term on the right-hand side of Eq. (2.7), after substituting Ut by the 

right-hand side of Eq. 113 becomes: 


e-27t 

L 


Re 



Utudx 


g-27t 

L 


k=l 


(A.l) 
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where 




Ji = —Re { j Uxxudx > , J2 = Re < i 
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J3 = Re ^ 7 j uudx > , J4 = Re 


— Re 


Jj — Re < 
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J {u - iaR){\u\‘^)xUudx 


-L 
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-L 

L 


{\uYu)xUdx > , 


Integration by parts and use of periodicity, implies that J5 = — J5, hence J5 = 0. By the same token, we 
have that j^j^{\u\^)x\u\^dx = 0, hence. 


\‘^u^udx 


Jq = /iRe / {\u\‘^)xuudx + /iRe 

J-L J-L 

/ L pL 

{\u\‘^)x\ufdx jil / \ufRe{uxu)dx 

-L J-L 

= fiRe J {\u\‘^):c\u\‘^dx + ^Re J {\u\‘^)x\u\‘^dx = 0, 

iu-iaR)j (|u7a;|u|^da;| = 0. 


Jj — Re < 


Similarly, we have that Ji = Re || = 0, and J 2 = Re |i = 0, as well. Thus, we 

recover that the only contribution t o (|a .i| comes from the terms J3 = 7 \u\^dx^ and J4 = \u\^dx. 

Inserting all the above relations to (lA.il), we observe that the latter is actually 




g-27t rL 

Re / utudx = ^—-— / \u\‘^dx-\-6 —-— / | 
J-L dj J_R L J_R 


—2jt pL 




(A.2) 


and consequently, Eq. (2.7) results in Eq. (2.8). 

Einally, we remark that following computations similar to the above, but for the (squared) L^-norm 
functional, we find that: 

pL pL pL 

(A.3) 


^ [ \u\‘^dx = 2j [ \u\^dx^25 f \u\^dx. 

dt J-L J-L J-L 


Thus, it is evident t hat th e factor exp(—27^) is the integrating factor for the linear part of the above equation. 
For 7 = d = 0, Eq. (A.3) is nothing but the conservation of the norm in the standard NLS model. 
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